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The train model which is a variant of the Burridge-Knopoff 
earthquake model is investigated for a velocity-strengthening 
friction law. It shows self-organized criticality with complex 
scaling exponents. That is, the probability density function 
of the avalanche strength is a power law times a log-periodic 
function. Exact results (scaling exponent: 3/2-|-27ri/ ln4) are 
found for a nonlocal cellular automaton which approximates 
the overdamped train model. Further the influence of random 
static friction is discussed. 

PACS numbers: 64.60.Lx, 05.70.Ln 

Ten years ago Bak, Tang, and Wiesenfeld showed that 
a weakly driven dissipative system with many metastable 
states can organize itself into a critical state in the sense 
of a second-order phase transition Because criti- 

cality is characterized by scale invariance they expect 
power-law behavior. For example the probability den- 
sity function of the strength S of the restructuring events 
(avalanches) should scale like where B is some pos- 

itive real number. 

In continuous scale invariance the scaling factor A can 
be arbitrary. This invariance is partially broken if A is 
restricted to a specific value A and its integer powers 
A" . This discrete scale invariance has a profound conse- 
quence It leads to complex scaling exponents B + iC, 
where C = 27r/lnA. More precisely: The scaling func- 
tion has the form S^^ f{ClnS), where / is a 27r-periodic 
function. Sornette and collaborators have shown that 
such scaling functions are common in many areas, e.g. 
fractals, deterministic chaos, dendritic growth and rup- 
ture §. 
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FIG. 1. The train model and the phenomenological friction 
law @. 

In this paper we present two models which exhibit self- 
organized criticality (SOC) with complex scaling expo- 
nents. The first model is the train model ^ which is a 
variant of the well-known Burridge-Knopoff (BK) earth- 
quake model Ql . In the literature the BK model and the 
train model are treated as examples of weakly driven dis- 



sipative systems with many metastable states. In both 
models power laws for the avalanche (earthquake) statis- 
tics have been found. In most of these studies an un- 
realistic purely velocity-weakening friction law is used. 
But any phenomenological friction law has to be velocity- 
strengthening for large velocities. Here we use a realistic 
friction law (see Fig. |^) which is a velocity independent 
Coulomb law for small velocities. For velocities larger 
than vq it is proportional to the velocity |^ . Our second 
model is a nonlocal cellular automaton which approxi- 
mates the train model in the overdamped limit. 

The train model is a finite chain of -I- 1 blocks on 
a rough surface. The blocks are coupled by springs and 
the interaction with the surface is described by a phe- 
nomenological dry-friction law F (see Fig. |l|): 

Mxj + F{xj) — K{xj+i — 2xj + Xj-i), j = 1, . . . , N, 

(1) 

where M is the mass of a block, Xj is the position of 
block j, and k is the stiffness of the springs. The system 
is driven by pulling block zero with a very small velocity 
V, i.e., xq = vt. The other end of the chain is free, i.e., 
xn+1 — xn. Our friction law F reads: 

r (-oo,Fs] ifi = 0; 
F{x) = < 'jvq if < i; < Wo; (2) 

[72; if i; > Wo- 

A resting block starts sliding if the sum of the spring 
forces is larger than F5. The friction law does not al- 
low backward motions because the static friction can 
take any negative value. The kinetic friction force is a 
monotonically increasing function of the velocity. We 
assume 7^0 < Fs otherwise the chain would not show 
avalanches. In the simulation we drive the system in- 
finitesimally slowly. This can be achieved in the follow- 
ing way. During an avalanche xq is held constant. After 
the avalanche, when all blocks are at rest (i.e., xj = 0), 
we set xo = {1 + e)Fs/K — X2 + 2xi, with e ^ 1. Thus 
the force on block number one is (1 -I- e)Fs and just ex- 
ceeds the threshold for sliding. Usually we have chosen 
£ = 10-4. 

Figure |^ shows the evolution of the model for two dif- 
ferent values of the damping constant 7. The avalanches 
always start at the pulling end. They propagate up to a 
certain block. This is the reason for the tree-like struc- 
tures seen in Fig. |[ We characterize the avalanches by 
two quantities which have a simple geometric meaning in 
Fig. ||: The length L and the strength S defined by 
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5' = ^|a;f°'^-a;j^°f°'-'=|, (3) 
i=i 

where x^''^°''° and are the position of the block j 

before and after the avalanche. The length L is the num- 
ber of blocks that are involved in the avalanche, i.e., 
^aftcr ^before ^ p^j, non-systcm-spanning avalanches 
(i.e., L < N) this quantity corresponds to the height 
of the branching points in the tree. In Fig. ^ the area 
below a branching point is just S. 





block position x. 



finite-size scaling ansatz, i.e., 

P{S) = S-^GiS/N"). (4) 

Note that B = a + 1 because the probability density 
function is the derivative of P. The inset of Fig. || shows 
the best approach to such an ansatz. The value of a 
was obtained from a fit of the averaged values of S for 
the system-spanning avalanches which are responsible for 
the last step in the cumulative density. This fit yields 
a = 2.34 ± 0.02. For the value of cr we have chosen 
a = 1 — 1/a because a sum rule (S) ^ N has to be 
fulfilled |]. 




FIG. 3. The cumulative density P{S). The parameters are 
7 = 0.1, vo = 0.01, M = Fs = K = l,a.ndN = 530, 750, 1060, 
1500, and 2121. The curves are shifted by 2 log A*' in order 
to separate them. The inset shows the results of finite-size 
scaling with a = 2.34 and a = 1 — 1/a = 0.573. 



FIG. 2. The positions of the blocks just before an avalanche 
starts. Several hundreds of avalanches are shown. The param- 
eters are iV = 50, vo = 0.01, M ^ Fs = n = 1. The initial 
values are Xj{0) = 0, for j = 0, . . . , A'^. 

We see a clear distinctive behavior between the un- 
derdamped case (7 = 0.1) and the overdamped case 
(7 = 4) The underdamped case is characterized 

by chaotic motion of the chain leading to an irregular 
sequence of avalanches even for nonrandom initial con- 
ditions 1^. In the overdamped case the motion is very 
regular. After the transient which is finished after the 
first system-spanning avalanche the same sequence of 
avalanches reappears periodically. In both cases SOC 
with discrete scale invariance occurs. But the details are 
different. 

First we investigate the underdamped case. Fig. ^ 
shows the cumulative density P{S) for five different val- 
ues of the system length N. The cumulative density P{S) 
is the probability to find an avalanche that is stronger 
than S. The cumulative densities for S show steps for 
large avalanches. Let us assume that P{S) fulfill a usual 



Although the finite-size scaling ansatz is not com- 
pletely satisfactory, one has the impression that for in- 
creasing N the scaling functions are approaching a saw- 
tooth function. In other words, the cumulative density 
has steps which becomes steeper and steeper for increas- 
ing N and larger S. The scaling function shows oscilla- 
tions which are periodic in the logarithm of the argument. 
This is a clear sign of discrete scale invariance. The nu- 
merically obtained scaling factor is A = 4.3 ± 0.4. Thus 
the complex exponent of the probability density function 
reads B + iC = a + 1 + 2ni/ \n\ ^ 1.573 + 4.3i. 

The selfsimilarity of the avalanche tree in Fig. ^(b) 
reflects the discrete scale invariance in the overdamped 
case. Between two system-spanning avalanches, 2" — 1 
avalanches occur. Here n = 1 -I- int (log A^/ log 2), where 
int (x) denotes the largest integer smaller than x. There 
are only n different avalanche lengths and sizes. The m- 
th type of avalanche occurs 2"^™ times. Its length is L = 
2m- 1 xhus, the scaling factor is two. The cumulative 
density P{L) is in a log-log plot a staircase with stairs of 
equal heights and widths |^ . A similar staircase is found 
for P{S). Here the scaling factor is A = 2^ = 4 because S 
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corresponds to an area in Fig. The critical exponents 
are a = 2 and a = 0.5. Thus B + iC = 3/2 + 27ri/ln4. 

The overdamped behavior of the train model can be 
mimicked by the following nonlocal cellular automaton. 
The state of each cell is given by the block positions Xj 
and a boolean variable Sj which is true when the block 
slides, otherwise it is false. The driving rule is the same 
as in the train model, i.e., xq = (1 + e)Fs / n — X2 + 2xi. 
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FIG. 4. The cumulative density P{S) for the deterministic 
cellular automaton for the overdamped train model with and 
without randomly chosen constant static forces Fs- The pa- 
rameters are (Fs) = k = 1 and A'' = 2121. The curves are 
shifted in order to distinguish them. The upper solid curve 
is the nonrandom case, i.e. AFs — 0. The middle dashed 
curve is one realization of the random static friction taken 
from a Gaussian distribution with AFs = 0.1. The lower 
solid curve is an average over 1000 realizations of the random 
static friction. The inset shows the scaling functions. 

The relaxation rules are the following: (i) First the 
forces fj = K{xj^i — 2xj + Xj+i) — Fs are calculated. 
The variable Sj is true if and only if /-,- > 0. (ii) In the 
second step the positions of all sliding blocks are updated 
simultaneously in the following way: 

Xj = Xj, + ^ ^(j - ji), for ji <j< j2, (5) 

J2 - Jl 

where ji and j2 are the left- nearest and right-nearest non- 
sliding blocks. That is, Sj^ = Sj^ — false and Sj = true, 
for jl < j < j2. The consequence of this rule is that 
Xj^i + Xj+i — 2xj — 0, for jl < j < j2- Note that 
So = s^v+i = false. If j2 = + 1 then — Xj^ in 
(|). This rule gives the result of a relaxation of sliding 
blocks governed by (|l|) assuming iiq = in the friction 
law (^) and immobile non-sliding blocks. The slightly 
curved lines in Fig ||(b) are the effect of vq ^ 0. (iii) 
In the third step the boolean variable Sj is recalculated: 
A sliding block still slides and nonsliding block starts to 
slide for the same reason as in the first step. That is, 

^now ^ ^oid ^ _ 2a;j + xj+i - Fs/k > 0), 

for j^l,...,N, (6) 



where V denotes the boolean operator for inclusive or. 
(iv) Repeat steps (ii) and (iii) until no new block starts 
to slide in step (iii). Note that rule (ii) is a nonlocal rule 
because ji and j2 can be arbitrary far from site j. This 
is in contrast to most other automata discussed in the 
field which have local relaxation rules. 

The result of applying these relaxation rules are the 
following. The first sliding block is always block number 
one. In the next cycle block number two starts to slide. 
This goes on until no new block starts to slide. Thus the 
avalanche length L is the smallest positive value which 
fulfills 

xo + (xL+i - xq)L/{L + 1) + XL+2 - 2xL+i < Fs/k. 

(7) 

After some straightforward calculations one finds that be- 
tween two system-spanning avalanches 2""™ avalanches 
occur. They are organized in the same binary tree as 
in the train model [see Fig. |^(b)]. The length and 
the strength of the m-th type are Lm — 2™"^ and 
Sm = (1 + 22™-1)Fs/(6k), respectively Therefore the 
cellular automaton has exactly the same scaling expo- 
nents and log-periodicity as the overdamped train model. 

The behavior of the automaton is very sensitive 
to nonuniformities in static friction Fs because the 
avalanches always stop when condition (|^) is just fulfilled. 
We introduce quenched randomness in two different ways 
leading to two different behaviors. 

In the first way random numbers Fsj from a Gaus- 
sian distribution are assigned to each block. After 
each system-spanning avalanche the same sequence of 
avalanches appears. They are also organized in a hierar- 
chical manner. The cumulative densities are still stairs 
(see Fig ^ but the heights and the widths of the steps 
are fluctuating. Averaging over many realizations of 
{Fsi,Fs2 ■ ■ ■ ,Fsn} leads to cumulative densities which 
still show log-periodic oscillations (see Fig. The os- 
cillation amplitude decreases with the noise level. The 
phase of the oscillation changes also but it still does not 
depend on N. 

The automaton shows a different behavior if we as- 
sign to each block a new random number Fsj after a 
slide. Figure || shows that the oscillations in the cu- 
mulative density vanish completely even for infinitesi- 
mally small noise level ||]. Otherwise the scaling ex- 
ponents are the same as for the nonrandom case. In the 
underdamped case the train model is less sensitive to 
this kind of quenched randomness. In simulations of (0) 
with AFs/{Fs) — 10^^ we stiU find oscillations but with 
smaller amplitudes. 

Why have these log-periodicities not be found for 
purely velocity- weakening friction laws |^? The main 
reason may be that for such laws any sliding motion 
of the chain or a part of the chain is unstable. For a 
velocity-strengthening friction force, sound waves with 
wavelength larger than L = AnV kM/^ are overdamped 
because their frequency becomes less than j/2M. Thus 
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avalanches involving more than L blocks show regular 
behavior after some possibly chaotic transients. This 
regular behavior is very well described by the nonlocal 
cellular automaton. In the overdamped case (i.e., L ^ 1) 
the regular motion successively amplifies the smallest in- 
trinsic length scale (i.e., L = 1) by the factor two. This 
basic mechanism is responsible for discrete scale invari- 
ance which would be destroyed by the intrinsic instabil- 
ities cause by a velocity-weakening friction force. It is 
unclear, why in the underdamped case (i.e., L'^ \) dis- 
crete scale invariance also occur. The mechanism has to 
be a different one because the motion of a group of less 
than L blocks is in general irregular. Two observations 
may help to explain this phenomenon: (i) The steps in 
the cumulative density of L do not occur at jDOwers of 
two but at aN/X'l, for n = 0, . . . , n, where Al ~ 2, a 
is an iV-independent number between 1/2 and 1, and 
h is clearly less than In (a7V/L)/ In 2. (ii) Although the 
large avalanches are organized in a binary tree, their se- 
quence of occurrence is different. Let us consider a big 
avalanche much bigger than L but not a system-spanning 
one. Looking for a considerably bigger avalanche in the 
past and in the future we find that in the overdamped 
case both time intervals are the same whereas in the un- 
derdamped case the time interval in the past is much 
smaller than the time interval in the future. In fact only 
a handful of very tiny avalanches occur in the time inter- 
val in the past. 




FIG. 5. The cumulative density P{S) of the cellular au- 
tomaton with randomly chosen static friction Fs which 
changes after each slide. The parameters are (Fs) — k — 1, 
AFc = 10"^ and = 530, 750, 1060, 1500, and 2121. The 
cumulative densities are shifted by 2 log N in order to separate 
them. The inset show the scaling function G. 

We have shown that SOC with complex scaling expo- 
nents occurs in the train model with a realistic velocity- 
strengthening friction law. We are confident that other 
depinning models may also show discrete scale invari- 
ance. The key ingredients are the following: (i) There 



should be a bistability between pinning and sliding for 
the same local force. This can be either achieved by in- 
ertia of the pinned objects or by age-dependent pinning 
(i.e., the pinning force which increases with the pinning 
time), (ii) The sliding dynamics should be nonchaotic 
with diffusion-like relaxation of long-wavelength excita- 
tions. In BK-like models a velocity-strengthening friction 
laws are a necessary condition for that, (iii) Quenched 
randomness in the pinning forces should be absent or 
weak. The first two properties are necessary to get non- 
local deterministic relaxations rules similar to our au- 
tomaton. These rules are responsible for discrete scale 
invariance because they either amplify successively an 
intrinsic (microscopic) length scale up to the system size 
by a constant factor or vice versa. Quenched randomness 
distort these process leading to fiuctuations of the phase 
of the log-periodicity which may be smeared out if they 
are too strong. 
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